#这部分代码是考虑双发和健康中没有8种pattern而作的分析。
setwd("E:\\5hmc_file\\2_5hmc_yjp_bam\\ASM\\bayes_pvalue_beta0")
group1=c("X2B_X1T","M8_M7","M6_M5","M2_M1","M48_M47","M50_M49","M28_M27","M30_M29","M26_M25","M35_M36","M18_M17","M20_M19","M22_M21","M40_M39")

i=1
fn1=paste0(group1[i],".bayes_p.txt")
file=read.table(fn1,head=T,sep = "\t")
file$unitID=paste(file$chrom,file$position,file$ref,file$var,sep=":")
file=file[file$normal_bayes_pvalue<0.05|file$tumor_bayes_pvalue<0.05,]
file=data.frame(unitID=file$unitID,normal_bayes_beta0=file$normal_bayes_beta0,normal_bayes_pvalue=file$normal_bayes_pvalue,tumor_bayes_beta0=file$tumor_bayes_beta0,tumor_bayes_pvalue=file$tumor_bayes_pvalue)
file$normal_group="nosig"
file[file$normal_bayes_beta0>0&file$normal_bayes_pvalue<0.05,]$normal_group="up"
file[file$normal_bayes_beta0<0&file$normal_bayes_pvalue<0.05,]$normal_group="down"
file$tumor_group="nosig"
file[file$tumor_bayes_beta0>0&file$tumor_bayes_pvalue<0.05,]$tumor_group="up"
file[file$tumor_bayes_beta0<0&file$tumor_bayes_pvalue<0.05,]$tumor_group="down"
file$pattern=paste0("normal_",file$normal_group,"-tumor_",file$tumor_group)
names(file)=c("unitID",paste(group1[i],names(file)[2:5],sep="_"),names(file)[6:8])
file$pattern2="other"###添加标签统计twins内的pattern
file[file$pattern=="normal_down-tumor_down",]$pattern2="normal_down_tumor_down"
file[file$pattern=="normal_down-tumor_nosig"|file$pattern=="normal_nosig-tumor_down",]$pattern2="one_is_down"
file[file$pattern=="normal_down-tumor_up"|file$pattern=="normal_up-tumor_down",]$pattern2="onedown_oneup"
file[file$pattern=="normal_nosig-tumor_up"|file$pattern=="normal_up-tumor_nosig",]$pattern2="one_is_up"
file[file$pattern=="normal_up-tumor_up",]$pattern2="normal_up_tumor_up"
rt=data.frame(table(file$pattern2))
names(rt)=c("Var1",group1[i])

###取并集
for(i in 2:length(group1)){
fn1=paste0(group1[i],".bayes_p.txt")
file=read.table(fn1,head=T,sep = "\t")
file$unitID=paste(file$chrom,file$position,file$ref,file$var,sep=":")
file=file[file$normal_bayes_pvalue<0.05|file$tumor_bayes_pvalue<0.05,]
file=data.frame(unitID=file$unitID,normal_bayes_beta0=file$normal_bayes_beta0,normal_bayes_pvalue=file$normal_bayes_pvalue,tumor_bayes_beta0=file$tumor_bayes_beta0,tumor_bayes_pvalue=file$tumor_bayes_pvalue)
file$normal_group="nosig"
file[file$normal_bayes_beta0>0&file$normal_bayes_pvalue<0.05,]$normal_group="up"
file[file$normal_bayes_beta0<0&file$normal_bayes_pvalue<0.05,]$normal_group="down"
file$tumor_group="nosig"
file[file$tumor_bayes_beta0>0&file$tumor_bayes_pvalue<0.05,]$tumor_group="up"
file[file$tumor_bayes_beta0<0&file$tumor_bayes_pvalue<0.05,]$tumor_group="down"
file$pattern=paste0("normal_",file$normal_group,"-tumor_",file$tumor_group)
names(file)=c("unitID",paste(group1[i],names(file)[2:5],sep="_"),names(file)[6:8])
file$pattern2="other"###添加标签统计twins内的pattern
file[file$pattern=="normal_down-tumor_down",]$pattern2="normal_down_tumor_down"
file[file$pattern=="normal_down-tumor_nosig"|file$pattern=="normal_nosig-tumor_down",]$pattern2="one_is_down"
file[file$pattern=="normal_down-tumor_up"|file$pattern=="normal_up-tumor_down",]$pattern2="onedown_oneup"
file[file$pattern=="normal_nosig-tumor_up"|file$pattern=="normal_up-tumor_nosig",]$pattern2="one_is_up"
file[file$pattern=="normal_up-tumor_up",]$pattern2="normal_up_tumor_up"
rtmp=data.frame(table(file$pattern2))
names(rtmp)=c("Var1",group1[i])
rt=merge(rt,rtmp,by="Var1")
}
write.csv(rt,"../20201104_at least one individuals is AShM/pvalue.sig.at.least.one.statis.5.pattern.csv",quote=F,row.names = F)

x=unlist(rt[1,2:7])
y=unlist(rt[1,12:15])
wilcox.test(x,y)